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Abstract 



We derive variationally the HFB equations for a trapped self-interacting Bose 
gas at finite temperature. In the Thomas-Fermi limit, we obtain simple expressions 
for the condensate, the non condensate and the anomalous densities. Their behavior 
in terms of the condensate fraction meets qualitatively the experimental data. In 
particular, the non condensate and the anomalous densities are peaked at the center 
of the trap and not at the edges as predicted by the self-consistent HFB calculations. 
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1 Introduction 



Since the discovery of the Bose-Einstein condensation, a great effort was devoted by 
researchers all around the world in order to understand and predict the condensate prop- 
erties. The main tools, beside the Monte-Carlo calculations [TJ, were the Bogoliubov[2], 
the Popov[3], the BeliaevjU [5] and the Hartree-Fock-Bogoliubov[6j El E] approxima- 
tions. These approximations all adopt some simplifying assumptions about the various 
order parameters such as the condensate density n c = |$| 2 , the non condensed density h 
and the anomalous density rh. 

In a previous paper [10J, we rely on the time- dependent variational principle of Balian 
and Veneroni[llj. which allows one to overcome some of the restrictions related to the 
various approximation schemes. We obtained a set of three coupled dynamical equations, 
which we called "Time-Dependent Hartree-Fock-Bogoliubov" (TDHFB) equations, gov- 
erning the evolution of the order parameter $, n and rh. They were shown to generalize 
in a consistent way the Gross-Pitaevskii equation[12j. 

In this paper, we focus on the static solutions of our TDHFB equations both at 
zero and finite temperature in the Thomas-Fermi (TF) regime. The interest is evident, 
since there remain many unanswered questions such as the general dependence of the 
density profiles on the temperature and on the interaction strength and the effect of the 
interactions on the critical temperature. More precisely, recent experiments are raising 
challenging questions about the structure of the thermal cloud and its backeffects on the 
condensate [13] and the role of the anomalous density. Two key questions are the observed 
compression effect of the thermal cloud on the condensate with increasing temperature 
and the fact that the thermal cloud and the anomalous average are peaked at the center 
of the trap. However, the full self-consistent HFB calculations [7J not only miss the first 
effect, but also predict that the non condensate and the anomalous densities are rather 
peaked at the edge of the condensate. One possible reason is the small atom numbers 
used in the above calculations. However, going to large atom numbers is clearly a quite 
formidable task in these methods. 

The goal of the present paper is to show that our equations reproduce correctly the 
compression effect and the shift of the maximum of the densities toward the borders of the 
trap when the number of atoms is large enough. Hence, the approximations we are actually 
using (mean field and Thomas- Fermi), when they are valid, retain the most important 
qualitative features without destroying the underlying physics. This may provide a simple 
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enough tool which can be considered as a starting point for a more elaborate treatment. 

After recalling the main steps that we have used to derive the variational TDHFB 
equations, we present in section 2 the static solutions and discuss their properties at 
zero temperature. At finite temperature, the equations are much more involved and re- 
quire a careful analysis. In the TF limit, we present a simple method which allows for a 
self-consistent determination of the various density profiles as well as some other static 
properties of the condensate such as the chemical potential and the condensate radius. In- 
deed, the TF approximation obviously provides simple enough analytical expressions since 
it amounts to neglecting the second order derivatives thus yielding algebraic equations 
instead of partial differential equations. This is the main advantage of our method which 
yields quite natural results without having to handle highly non-linear self-consistent 
equations. 

In section 3, we present the results of our calculations. We plot first the condensate 
radius and the central density as functions of the condensate fraction and note in particular 
the compression effect of the condensate due to the thermal cloud. Moreover, we discuss 
the TF profile obtained for the condensate density even at low condensate fraction. The 
non condensate density is also plotted for a wide range of condensate fraction and shows 
a good qualitative agreement with recent experiments. Finally, the anomalous density, 
although not yet measured experimentally, is shown to behave in a quite intuitive way. 

Some concluding remarks and improvements of our method are given at the end of 
the paper. 

2 The Variational HFB Equations and their Thomas- 
Fermi Limit 

The general TDHFB equations were derived in ref . [TO] for a grand canonical Hamiltonian 
of trapped bosons with mass m and quartic self-interactions (with coupling constant g): 

o(r) + - / a+(r)a + (r)a(r)a(r). (2.1) 

2 Jr 

The quantity V ext (r) is the trapping potential and fi is the chemical potential. These 
equations read: 



H 



h 2 

-—A + V ext (r)-fi 
2m 
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ih ® = ( _ fi A + y ext -v + gn c + 2gn) $ + #m$*, 

ihn = g (m*$ 2 - m$* 2 ) , (2.2) 
ihfh = ±(2Vn + 1)$ 2 + 4 (-J^A + Kxt - /i + 2#n + &(2Vn + 1)) m, 

where V is the volume of the gas. In Eqs. ( 12.21) . $ is the order parameter, n c the conden- 
sate density (n c = |$| 2 ), n the non-condensed density (or thermal cloud) and rh is the 
anomalous density. The quantity n = n c + n is the total density. 

The TDHFB equations with a general Hamiltonian H were derived in [T4]. The 
properties discussed here and in [10] were established there in their general forms. These 
equations were obtained using the Balian-Veneroni variational principle [HJ, with a gaus- 
sian trial density operator (that is, an exponential operator of a quadratic form) in the 
creation and annihilation operators. The result was a set of coupled evolution equations 
for the expectation values (a), (a + a) — (a + )(a) and (aa) — \ (a}\ 2 . When one identifies 
these quantities respectively with the order parameter the non-condensed density h 
and the anomalous density rh, and when one restricts H to the class (12. ip . the equations 
( Q follow. 

The TDHFB equations couple in a consistent and closed way the three densities. 
They should in principle yield the general time, space and temperature dependence of 
the various densities. Furthermore, they obviously constitute a natural extension of the 
Gross-Pitaevskii equation[12]. They are not only energy and number conserving, but 
also satisfy the Hugenholtz-Pines theorem (see below) which leads to a gapless excitation 
spectrum in the uniform limit. The two last equations in ( 12. 2p are not totally independent 
since n and rh are related by the "unitarity" relation[10j: 

I = (l + 2Vnf -(2V\m\)\ (2.3) 

where the Heisenberg parameter / (which is always > 1) is a measure of the temperature, 
the lower limit (J = 1) being the zero temperature case. For instance, for a thermal 
distribution at equilibrium, / writes as / = coth 2 (huo/2kBT), where ujq is the average 
frequency of the trapping field [10] 0- We therefore see that upon replacing n by its 
expression given in fl 2 . 3 j) . the temperature appears explicitly in the equations. 

It is to be mentioned that the TDHFB equations have also been derived by other 
groups using different variational f ormulat ions 1 1 5 1 [16]. For instance, the authors of the 
first reference have used a gaussian trial wave functional and have obtained a set of 

1 In fact, one can show that for a system of energy E, \fl =1 + 2 Jb(E), where /s is the Bose-Einstcin 
distribution. 
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equations very similar to ours. However, their discussion is purely formal and does not 
address for instance the question of validity of their approximation. 

The static solutions, which are the object of our study in this work, are obtained 
by setting to zero the right hand sides of ( 12. 2p . At zero temperature, the standard TF 
limit pE] amounts to neglecting the kinetic (or A) term in the Gross-Pitaevskii equation. 
This is particularly satisfied for strong interacting regimes or large atom numbers. At 
finite temperature and below the transition, since there are two phases (condensed and 
non condensed) which coexist, one has to provide a complementary recipe for what we 
shall call the finite temperature TF limit. First, neglecting the kinetic energy of the 
condensate remains a justifiable approximation since the atoms are slowed down in order 
to obtain condensation. On the other hand, rh is believed to be an extremely small 
and slowly varying function whatever the temperature is (recall that it describes the 
correlations between the condensed and non-condensed phases). Hence, one may in a first 
approximation safely neglect Am. Heuristically, one may argue that, since the equations 
for n c and rh contain almost comparable operators, h Q and h + g(n c + n/2), where ho is 
the self-consistent mean field hamiltonian h = Kxt( r ) — A 4 + 9 n c + 2gn, the TF condition 
ho » T (T being the kinetic operator), if fulfilled for n c should also be satisfied for rh. 
For this approximation to be consistent, n c and rh should vary on the same characteristic 
length, which is indeed the case as we will show later. 

Before proceeding further, it is important to notice at this point that a kinetic-like 
term of the thermal cloud does not appear explicitly in the equations but is rather hidden 
in the third equation of ( 12.21) . Indeed, the kinetic term of the thermal cloud is related to 
the second derivative of the anomalous density. Differentiating (12.31) yields a relation of 
the form: 

Ah ~ (V|m|) 2 - (V|n|) 2 + |m|A|m|, (2.4) 

which shows in particular that neglecting Am does not necessarily mean neglecting An. 
That is precisely the recipe that we shall adopt below. 

With this finite temperature prescription, the static equations now write 
(Kxt(r) - u + gn c + 2gn) $ + grh<5>* = 0, 

m*$ 2 - m$* 2 = 0, (2.5) 
(Vcxt(r) - u + 2gn + ^{2Vh + 1)) rh + ^(2Vn + 1)$ 2 = 0, 

These equations are naturally gapless and satisfy the Hugenholtz-Pines theoremjBj. In- 
deed, owing to the second equation in (12.51) . one can easily show that at zero momentum, 
the relation ii — g(n + h — \rh\) is clearly satisfied without adding further assumptions, 
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as is usually performed [6]. 

In order to solve these equations, we may distinguish two different situations. The 
first one is for T = 0. When all the atoms are condensed, n = fh = 0, and n c equals the 
total density n of the gas. Omitting the trivial solution with n c = 0, one may take into 
account just the first equation in (12.51) . since we consider a gas without a quantum cloud. 
Indeed, within the present set of equations, it is an approximation (although justifiable) to 
ignore the quantum depletion at T = 0. The last two equations in (12. 5p become therefore 
meaningless, and we are left with a simple expression for the condensate density 

n c (r) = -£(r) = -(//-V ext (r)). (2.6) 
9 

Upon defining the oscillator length a = (h/mcuo) 1 ^ 2 and the s-wave scattering length a = 
mg/ATrh 2 , we obtain for a spherical trapping potential V ext (r) = ^muj 2 r 2 , the condensate 
radius R and the reduced chemical potential z/ = \ij \hu Q for a gas of N bosons as 

R { a\ 1 / 5 

W-) , (2.7) 



ao V ao 



a \ 2 / 5 

m=(l5JV— ) . (2.8) 



The preceding expressions show that the spreading of the condensate depends essentially 
on the balance between the self-interactions and the trapping potential. These results 
have also been obtained by many other authors, see e.g. [5j l9l fTU]. 

When < T < T BEC , we have of course n c ^ and h ^ 0. Let us introduce the 
parametrization 2Vn + 1 = vT cosher, 2V|m| = A/Tsinha, which automatically endows 
the relation (12.31) . Then, from the third equation in (12.51) . one obtains a simple equation 
for X = e a : 

3X 4 - AX 2 + 1 + ^p- (X 2 - 3) X = 0, (2.9) 



yfl 

from which one extracts n and \m\ as functions of n c . Next, one uses these expressions in 
the first equation (12.51) to get the condensate density 

n c (r) = -£(r) - i ( - ^ ' ^ ' V7 - l) . (2.10) 

What is remarkable is that the sole acceptable solution of equation (12.91) is a bounded 
function of rj = Vn c /yl. It is represented on figure 1. 
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Fig.l: Solution of Eq. (l2.9p vs. the dimensionless condensate density. 



Due to this behavior, one can easily show that the quantity 



which appears 

in (12.1 Op is almost independent of n c and becomes rapidly close to unity. Hence, one may 
safely approximate (I2.10p by 



n c {r) ~ -f (r) - — (v 7 / - l) . 



(2.11) 



In fact, one can check that the relative error between the two expressions (12.101) and 
(12.111) is less than 1%. Finally, since vT does not depend on space, the result ( 12.111) 
shows that the finite temperature correction to the Thomas-Fermi profile (12.61) is simply 
a space-independent (but temperature dependent) shift. The condensate density finally 
writes in the suggestive form 



n c (r) 



Kxt(fl)~Kxt(r) 
9 



(2.12) 



which is formally the zero temperature TF profile. It is then easy to show that the 
condensate radius R takes also a simple form 



— = I r,\ — 



1/5 



(2.13) 



do V CLq , 

but now, it is the number of condensed atoms N c which is involved and not the total 
number of atoms. Hence, our finite temperature prescription for the TF approximation 
provides a natural extension of the zero temperature expressions, since the Thomas-Fermi 
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parameter is now N c a/ao instead of Nd/a®. The same conclusion may be drawn for the 
chemical potential which now writes 

v= 15iV c — +6tt— (V7-1 15iV c — . (2.14) 
V a J a \ a J 

In order to apprehend better these results, let us compute the remaining unknown 
quantities, such as the non condensed and the anomalous densities. To this end, and in 
order to obtain tractable expressions, we find it more convenient to use the simple fit 

X = , (2,5) 

77 + 2/3 

instead of the full analytical solution of equation (12.91) which reproduces correctly the 
solution X plotted in figure 1 with a residual error less than 0.1%. Upon rewriting 
equation (12.121) in the form 77 = 770 (1 — x 2 ), with an obvious definition of 770 , we readily 

get 

1 I VI fv^o(l-f) + 2/3 .od-^) + 2/3 \ 1 
KJ 2V\ 2 V 770(1 - a: 2 ) + 2/3 ^ (1 - x 2 ) + 2/3/ J 1 ' 

, , / n = J_ VI / v / 3^o(l-a: 2 )+2/3 _ 77 (1 - x 2 ) + 2/3 \ 

|m|W 21/2 Uo(l-, 2 ) + 2/3 V377o(l-x 2 )+ 2/3 J' 1 J 

where x is the reduced radial distance x = r/R. 

In order to obtain more quantitative results, one must determine N c by using the 
normalization condition. We get easily the relation 

1 + 2N = 2N C + Vlg(s), (2.18) 

where 



2 . r ( 3 , , 1 

—j= + (V3 — l)s < 1 Vs + larctanh 

V 3 [ V s + 1 



q(s) = —= + (V3 — l)s<(l V s + larctanh . H — J —= + larctanh . 

(2.19) 

with s = 4v7/15iV c . But since the function g(s) satisfies 1 < g(s) < 2/v3, the equation 
(12.181) is approximately solved to yield, to a very good accuracy, the simple result 

N c ~N-^^±. (2.20) 

All the unknown quantities may now be determined in terms of N and y/l alone. The 
corresponding results will be discussed in the next section. 
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3 Results and Discussions 



First of all, the condensate radius (I2.13P may be written as 




) 



1/5 



(3.1) 



where -Rtf is the zero temperature result given by equation (12.71) . Figure 2 represents 
the condensate radius (in units of -Rtf) as a function of the condensate fraction and 
we notice in particular the compression of the condensate when reducing N c /N (that is 
increasing the temperature). This effect is by now a well established experimental result 
[T3] and is attributed to the thermal cloud. The same effect of compression is observed 
on figure 3 for the central density n c (r = 0) but it is more pronounced due to the power 
law of 2/5 (see 12.12"]) instead of 1/5 for the condensate radius. To be more precise, let 
us choose generic values for the number of atoms and the interaction strength (N = 10 5 
and a/a = 0.5 10~ 3 ) and plot the various densities (in units of the oscillator volume a%) 
versus the radial distance (in units of .Rtf = 3.758a ) for a condensate fraction ranging 
from 5% up to 60%. 
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Fig. 2: Condensate radius vs. the condensate fraction. 
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N c /N 



Fig. 3: Central density vs. the condensate fraction. 

The figure 4 shows typical Thomas- Fermi profiles for the condensate density, even 
for low condensate fraction. This is of course what one may expect on general grounds 
in the TF regime. Moreover, the effect of compression of the condensate is also clearly 
visible here. 



1000 -, 




Fig. 4: Condensate density (in units of a 3 ) vs. the radial distance (in units of Rtf, see 

text) for N = 10 5 and a/a = 0.5KT 3 . 
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The non condensate density (from which we have subtracted a constant n(R) for clarity) 
is plotted on figure 5 with the same units as before. 

NJN=5% 




Fig. 5: Non condensate density (in units of a 3 ) vs. the radial distance (in units of -Rtf) 

for N = 10 5 and a/a = 0.51(T 3 . 

We notice that when increasing the condensate fraction, the thermal cloud tends to 
spread and flatten. This behavior is nicely comparable to the figure 4 of Ref . [13] (we have 
indeed deliberately chosen the same values for the condensate fraction. See however the 
remarks in sect. V of [13^). Furthermore, the non condensate density is peaked at the 
center of the trap for the whole range of the condensate fraction that we have considered. 
This is also confirmed by the experimental results of Gerbier et al. On the other hand, the 
full self-consistent HFB calculations [7] predict that the non condensate density is sharply 
peaked at the edge of the condensate, having therefore a "hole" at the center of the trap. 
It seems however that this structure would shift toward the center with increasing atom 
number, but this behavior has not yet been checked [TU]. 

Nevertheless, we do understand in our case that neglecting the kinetic-like terms 
entails overestimating the (effective) repulsive self-interactions which tend to make the 
profiles more uniform. This overall effect is also observed on the anomalous density (fig. 
6). Even if this quantity has not yet been measured experimentally, it is interesting to 
notice that our calculations predict a very simple and yet intuitive behavior as well which 
has been confirmed elsewhere [7]. 

Although, in order to compare exactly with the experimental results, one must include an overall 
scale factor due to the finite ballistic expansion time. 
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On the other hand, the thermal cloud takes on a (small but) finite value for r > 
R. Even if this behavior is less intuitive, it is not very surprising since we do know 
that neglecting the second order derivatives amounts to making a cut of the densities 
at the boundaries. It is indeed a limitation of the TF approximation as a whole at the 
boundaries [18]. The tail should be reproduced when one reinjects the second derivatives 
in the equations. In fact, it has been shown in the context of the full HFB numerical 
calculation [7] that this is indeed the case. 



m(r) 
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Fig. 6: Anomalous density (in units of a 3 ) vs. the radial distance (in units of -Rtf) for 

N = 10 5 and a/a = 0.5HT 3 . 

Finally, we observe on the figures 4, 5 and 6 that n c , \m\ and n vary on the same length 
scale (R) which justifies a posteriori our previous assumption. 



4 Concluding Remarks 

We present in this paper a finite temperature analysis of the HFB equations in the 
Thomas-Fermi limit for a gas of bosons in a harmonic trap. At zero temperature, we 
obtain familiar expressions for the chemical potential and the condensate radius. The 
standard Thomas-Fermi profile for the condensate density is also recovered. 

At finite temperatures and below the transition, since there are two phases, one 
should provide a prescription for the TF limit. We propose such a recipe (maybe the sim- 
plest) which consists in neglecting the second order derivatives of the condensate density 
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and the anomalous density. The underlying idea is that, although the anomalous density 
is necessary for the coherence of the equations, it is believed to be a very small and a very 
smooth quantity. We therefore obtain analytical expressions for the condensate density, 
the condensate radius, the chemical potential and the central density as functions of the 
condensate fraction (or the temperature). Our expressions appear as natural extensions 
of the zero temperature TF limit. 

Most importantly, we derive quite simple expressions for the non condensate density 
and for the anomalous density, which we plot as functions of the condensate fraction and 
draw many conclusions. First of all, the condensate profile is almost of the TF shape of 
which the spatial extension and the heights are controlled by N c /N . Furthermore, the 
compression of the condensate by the thermal cloud with increasing temperature is clearly 
visible. On the other hand, the non condensate density profile is qualitatively consistent 
with the temperature dependence observed in recent experiments. In particular, the non 
condensate density is not peaked at the edge of the condensate (as predicted by the 
HFB calculations) but has a broad maximum at the center of the trap. The thermal 
cloud tends to spread and flatten with increasing temperature. The calculated anomalous 
density, although not yet observed experimentally, shows also a very intuitive behavior; 
it has a broad maximum at the center of the trap and vanishes at the boundaries. The 
tendency to spreading and flattening with increasing temperature is also observed here. 

At the borders of the trap (where the condensate density vanishes) and for a given 
temperature, the non condensate density takes on a finite value which is a quite abrupt 
behavior. Although this meets the fact that the thermal cloud is actually surrounding 
the condensate, it is by no means conclusive. But this is a known shortcoming of the 
TF approximation as a whole since it breaks down at the boundaries of the condensate. 
Indeed, reinjecting the second derivatives of the densities will entail a more physical 
behavior. 

These results are quite satisfactory considering the drastic simplifications that occur, 
as compared for instance to the full self-consistent HFB calculations. Indeed, even if our 
overall density profiles are smoother and broader than what they should be, we know 
precisely the reasons and how to incorporate gradually what is missing. 

It is important to emphasize that the self-consistent HFB equations obviously cannot 
be easily handled for high atom numbers. Hence, our equations may provide a simple 
tool to gain new insights in situations where these techniques become computationally 
consuming. One may make contact with the HFB calculations by taking into account 
the kinetic terms. That is precisely what we are actually performing in the semi-classical 
limit. The preliminary results show a qualitative and quantitative improvement of the 
predictions. More details are postponed to a forthcoming paper [20J. 



13 



We are grateful to P. Schuck and Y. Castin for fruitful discussions and a careful 
reading of the manuscript. Special thanks to the members of the Groupe de Physique 
Theorique, IPN-Orsay-France, where part of this work has been done. 



14 



References 



[1] W. Krauth, Phys. Rev. Lett. 77, (1996), 3695; A. B. Kuklov, N. Chencinski, A. M. 
Levine, W. M. Schreiber and J. L. Birman, Phys. Rev. A55, (1997), 488; N. V. 
Prokof'ev, B. V. Svistunov and I. S. Tupitsyn, Sov. Phys. -JETP87, 310 (1998). 

[2] N. Bogoliubov, J. Phys. USSR 11 (1947), 23; A. L. Fetter and J. D. Walecka, "Quan- 
tum Theory of Many-Particle Systems", McGraw-Hill, NY, 1971; L. Pitaevskii and 
S. Stringari, " Bose-Einstein Condensation", International Series of Monographs on 
Physics, Oxford Science Publications, Clarendon Press, Oxford, 2003; C. J. Pethick 
and H. Smith, "Bose-Einstein Condensation in Dilute Gases", Cambridge University 
Press, Cambridge, UK, 2002. 

[3] V. N. Popov, Sov. Phys. JETP 20 (1965), 1185; "Functional Integrals and Collective 
Excitations", Cambridge Univ. Press, Cambridge, 1987; D. A. W. Hutchinson and 
E. Zaremba, Phys. Rev. A57 (1998), 1280; S. A. Gardiner, Phys. Rev. A56 (1997), 
1414; S. A. Morgan, J. Phys. B33 (2000), 3847. 

[4] S. T. Beliaev, Soviet. Phys. JETP 7 (1958), 289. 

[5] A. Griffin and H. Shi, Phys. Rep. 304 (1998), 1. 

[6] A. Griffin, Phys. Rev. B53, (1996), 9341; N. P. Proukakis and K. Burnett, J. Res. 
Natl. Inst. Stand. Technol. 101, (1996), 457; D. A. W. Hutchinson, E. Zaremba and 
A. Griffin, Phys. Rev. Lett. 78, (1997), 1842; D. S. Jin, M. R. Matthews, J. R. 
Ensher, C. E. Wieman and E. A. Cornell, Phys. Rev. Lett. 78, (1997), 764; V. M. 
Perez-Garcia, H. Michinel, J. I. Cirac, M. Lewenstein and P. Zoller, Phys. Rev. A56, 
(1997), 1424; D. A. W. Hutchinson, R. J. Dodd and K. Burnett, Phys. Rev. Lett. 81 
(11), (1998), 2298; T. Kita, J. Phys. Soc. Jpn 75, 044603 (2006); V. I. Yukalov and 
H. Kleinert, Phys. Rev. A75, 063612 (2006). 

[7] D. A. W. Hutchinson, R. J. Dodd, K. Burnett, S. A. Morgan, M. Rusch, E. Zaremba, 
N.P Proukakis, M. Edwards and C. W. Clark, J. Phys. B33, (2000), 3825. 

[8] A. K. Kerman and P. Tommasini, Phys. Rev. B56, (1997), 14733; Ann. of Phys. 
(NY.) 260 (1997), 250. 

[9] F. Dalfovo, S. Giorgini, L. P. Pitaevskii and S. Stringari, Rev. Mod. Phys. 71, (1999), 
463; J. Javanainen, Phys. Rev. A54 (1996), 3722; A. Minguzzi and M. P. Tosi, J. 
Phys.: Condens. Matter 9 (1997), 10211; S. Giorgini, Phys. Rev. A57 (1998), 2949. 

[10] M. Benarous, Ann. of Phys. (N.Y.) 320 (2005), 226. 



15 



[11] R. Balian and M. Veneroni, Ann. of Phys. (N.Y.) 187 (1988), 29; Ann. of Phys. 
(N.Y.) 195 (1989), 324. 

[12] E. P. Gross, Nuovo Cimento 20 (1961), 454; L. Pitaevskii, Soviet Phys. JETP 13 
(1961), 451. 

[13] F. Gerbier, J. H. Thywissen, S. Richard, M. Hugbart, P.Bouyer and A. Aspect, Phys. 
Rev. A70, (2004), 013607. 

[14] M. Benarous and H. Flocard, Ann. of Phys. (N.Y.) 273 (1999), 242. 

[15] V. Chernyak, S. Choi and S. Mukamel, Phys. Rev. A67 (2003), 053604. 

[16] L. Pricoupenko, Phys. Rev. A70 (2004), 013601. 

[17] G. Baym and C. J. Pethick, Phys. Rev. Lett. 76 (1996), 6; A. L. Fetter, J. Low 
Temp. Phys. 106 (1997), 643; E. Timmermans and P. Tommasini, Phys. Rev. A55 
(1997), 3645; S. Stenholm, Phys. Rev. A57 (1998), 2942; P. Schuck and X. Vinas, 
Phys. Rev. A61 (2000), 43603. 

[18] F. Dalfovo, L. P. Pitaevskii and S. Stringari, Phys. Rev. A54 (1996), 4213; A. L. 
Fetter and D. L. Feder, Phys. Rev. A58 (1998), 3185; A. Muhoz Mateo and V. 
Delgado, Phys. Rev. A74 (2006), 065602. 

[19] D. A. W. Hutchinson, Private communication. 

[20] M. Benarous and P. Schuck, In progress. 



16 



